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Abstract 

Terrestrial Gamma-ray Flashes (TGFs) are short intense flashes of gamma rays associated with 
lightning activity in thunderstorms. Using Monte Garlo simulations of the Relativistic Runaway 
Electron Avalanche (RREA) process, theoretical predictions for the temporal and spectral evo¬ 
lution of TGEs are compared to observations made with the Gamma-ray Burst Monitor (GBM) 
on board the Fermi Gamma-ray Space Telescope. Assuming a single source altitude of 15 km, 
a comparison of simulations to data is performed for a range of empirically chosen source elec¬ 
tron variation timescales. The data exhibit a clear softening with increased source distance, in 
qualitative agreement with theoretical predictions. The simulated spectra follow this trend in the 
data, but tend to underestimate the observed hardness. Such a discrepancy may imply that the 
basic RREA model is not sufficient. Alternatively, a TGF beam that is tilted with respect to the 
zenith could produce an evolution with source distance that is compatible with the data. Based 
on these results, we propose that the source electron distributions of TGEs observed by GBM vary 
on timescales of at least tens of microseconds, with an upper limit of ~100 ^s. 
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I. INTRODUCTION 


Terrestrial Gamma-ray Flashes (TGFs) are short intense flashes of gamma rays associated 
with lightning activity in thunderstorms which were discovered serendipitously in 1994 by 
the Burst And Transient Source Experiment (BATSE) pQ. TGFs are characterised by short 
timescales (< 1ms) and hard spectra which can extend up to tens of MeV [2]. Since their 
discovery, TGFs have been extensively studied by the Reuven Ramaty High Energy Solar 
Spectroscopic Imager (RHESSI) |3], Astrorivelatore Gamma a Immagini LEggero (AGILE) 
[1], and the Gamma-ray Burst Monitor (GBM) and Large Area Telescope (LAT) on board 
the Fermi Gamma-ray Space Telescope The exact emission mechanism of TGFs 

is unknown, but the leading theoretical models involve the Relativistic Runaway Electron 
Avalanche (RREA) process, whereby electrons are accelerated to high energies in electric 
fields [7]. As they propagate through the atmosphere, these electrons emit gamma rays via 
bremsstrahlung. The spectral and temporal properties of many averaged TGFs have been 
compared to RREA simulations and found to be broadly consistent for RHESSI observations, 
e.g. Smith et al. [2], Grefenstette et al. [8], Marisaldi et al. |9]. However, observations of a 
power law extending up to 100 MeV by AGILE have challenged this view, as such a spectral 
shape is inconsistent with standard RREA models 

The high count rate (on the order of hundreds of kHz) and low statistics associated with 
TGFs greatly complicates their analysis. The large effective area of the BATSE detectors 
allowed the study of TGFs on an individual basis, e.g. Feng et al. HI], Nemiroff et al. 
na, Gstgaard et al. PI- However, these observations were later found to have been heavily 
modified by instrumental dead time mm- In general, RHESSI does not collect enough 
counts per TGF to study them on an individual basis. Gonsequently, analysis of this data 
has concentrated on stacking many TGFs and studying the average behaviour, e.g. Smith 
et al. [2], Grefenstette et al. [3]. The stacking of RHESSI events is necessary but unfortunate, 
as it combines many TGFs with differing orientations and source - detector geometries. As 
the distance between the source of the TGF and the observer is increased, the entire TGF 
is expected to soften and be temporally extended, as a greater proportion of the collected 
counts have undergone Gompton scattering [I31[I51[I6] due to the greater integrated depth 
of atmosphere traversed. 

TGFs typically consist of an individual pulse, but can be composed of multiple emission 
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episodes. The time profile of a pulse can be either symmetric (comparable rise and fall time) 
or asymmetric, with a faster rise time than fall time HT]. This temporal asymmetry and 
the tendency to soften in time was hrst noted in BATSE TGFs [12]. This softening was 
quantihed as the temporal lag between the peak of soft (25-110 keV) and hard (> 110 keV) 
counts in 15 BATSE TGFs m- The lags ranged from ~ 70-370 ps, with an average value 
of ~100 fis. A detailed analysis of Monte Carlo simulations showed that these lags could be 
explained as being solely due to Compton scattering |T3|. A study of 431 combined dead 
time corrected RHESSI TGFs found an average delay between the arrival times of the soft 
(<300 keV) and hard counts to be 28 ±3 ps [S]. The discrepancy between the two values is 
due to the effects of dead time which were not accounted for in the initial studies of BATSE 
data. 

Although not as sensitive as BATSE, GBM still accumulates enough counts per TGF 
that they can be studied on an individual basis. The implications of dead time on GBM 
observations are also less severe than BATSE. The absolute timing accuracy of GBM allows 
correlations of TGFs detected with individual lightning strokes, e.g. Connaughton et al. 
[I8] . For the subset of TGFs detected by GBM with associated radio detections, the source 
location and hence the distance and relative orientation of the spacecraft to the source 
can be determined. This allows a more detailed study of the temporal evolution of TGFs 
as a function of source to satellite distance than has previously been possible. Using the 
Runaway Electron Avalanche Model (REAM), a Monte Carlo code developed by J. R. Dwyer 
at Florida Institute of Technology to simulate the RREA model m 120] . we obtain the 
predicted temporal evolution as a function of source to satellite distance. However, these 
predictions cannot be directly compared to the observations of GBM, as the counts observed 
by an instrument have been distorted by the detection process. Additional distorting effects 
can arise due to instrumental effects such as dead time and pulse pile-up. 

In order to compare the model predictions to the observations, the predictions must be 
folded through the Detector Response Matrices (DRMs) and passed through a dead time 
hlter. A GBM DRM is a lookup table which maps incident photons to detector counts. It 
is not diagonal, as incident photons do not necessarily have to deposit all their energy in 
the detector. The DRMs and dead time distort the incident photon distribution, with the 
result that varying source models may appear similar. This degeneracy is unfortunate, but 
unavoidable with the current generation of instruments. 
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II. OBSERVATIONS AND DATA 


The Fermi Gamma-ray Space Telescope consists of two instruments, the Large Area 
Telescope (LAT) |2T] and the Gamma-ray Burst Monitor (GBM) |22]. GBM consists of 14 
individual scintillation detectors, 12 sodium iodide (Nal) with an effective energy range of 
10-1000 keV, and two bismuth germanate (BGO) crystals with an effective energy range of 
0.2-40MeV. The large effective area of the BGO detectors (~ 160 cm^ per detector |23j ) 
is particularly important as it facilitates the accumulation of sufficient statistics such that 
TGFs can be studied on an individual basis. 

The relative timing resolution of 2 /rs and absolute timing accuracy of several microsec¬ 
onds, allows the correlation of gamma rays with individual lightning strokes. This is vital, 
as it allows the source location, and therefore the orientation and off-axis distance of the 
spacecraft relative to the source to be calculated. 

In the period from the launch of Fermi to March 16*'^ 2013, GBM has detected ~ 1330 
TGFs. This includes TGFs which triggered GBM and those which have been found in 
an off-line search. A preliminary correlation study has confirmed radio matches with the 
World Wide Lightning Location Network (WWLLN) and/or the EARTH Networks Total 
Lightning Network (ENTLN) for 287 of these TGFs. For each TGF in the sample, the radio 
location was used to determine the distance and orientation of Fermi to the source. Using 
this information, DRMs were generated for each individual TGF. The data used in this work 
are GBM Time-Tagged Events (TTE), which have a relative timing resolution of 2.0 /rs and 
128 pseudo-logarithmically spaced energy channels. 


III. SIMULATIONS 

Gomparing RREA simulations to individual TGFs is an extremely complex task, as the 
exact electric held {E) orientation and beam geometry at the source is unknown. However, a 
statistical study in which we assume an altitude and beam geometry can be used to study the 
effect of increasing source offset from the satellite on the simulations. These simulations can 
be folded through the DRMs and dead time hltered to create ‘synthetic’ TGFs. These can 
then be analysed in exactly the same fashion as the data. In this way, a direct comparison 
can be made between the theoretical predictions and the observed trends in the data. 
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A. Simulating RREA time profiles 


To generate the predicted temporal and spectral distributions, we use the REAM Monte 
Carlo simulation code with an ambient DC E held value of 400 kV/m. An instantaneous 
electron distribution of the form 

/(E) ocexp(-E/7.3MeV) (1) 

is created at the source with 10^ electrons. In the presence of the E held, these electrons 
undergo RREA multiplication. The electron interactions in the simulation include ionisa¬ 
tion, atomic excitation and Mpller scattering. Elastic scattering is fully modelled using a 
shielded Coulomb potential and includes bremsstrahlung production of X-rays and gamma 
rays. Compton scattering is also fully modelled and allows for the production of secondary 
electrons. After hve avalanche lengths (A ~ 50m), the resulting distribution of ~ 1.6 x 10® 
photons is propagated from the source to an altitude of 100 km. At this point, the atmo¬ 
spheric density is sufficiently low that the photons can be simply translated to 565 km, the 
altitude of GBM. These electrons travel along held lines and can also create more photons 
via bremsstrahlung. 

Once propagated to the spacecraft altitude, the photon distribution is integrated into 
concentric annuli of diameter 50 km. This ensures that there are sufficient statistics in each 
ring while also allowing the ehect of an increased source distance to be studied. We use 
a source altitude of 15 km and a wide beam geometry, as this is believed to be typical 
for TGFs (e.g. Smith et al. [2], Grefenstette et al. [H], Ostgaard et al. [13], Garlson et al. 
[24] . Hazelton et al. [23], Gjesteland et al. [23] )• For this geometry, the photons are spatially 
spread into a 45 degree isotropic cone (half-angle), simulating a diverging electric held at the 
source. Figure shows the evolution of the temporal and spectral properties as a function 
of source distance. In this Figure, the expected spectral softening and temporal elongation 
with increased distance is clearly evident. 

B. Synthetic TGFs 

To generate synthetic TGFs with RREA time prohles and spectra, we hrst consider a 
source directly below the satellite. This source is moved away from the satellite in 50 km 
steps. For each step, 1000 synthetic TGFs are generated by randomly drawing a huence 
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of photons for each detector from the predicted temporal and spectral distribntions for 
that offset. These are then folded throngh the DRMs and a dead time hlter is applied. 
Poisson noise is added to each detector to simnlate a backgronnd component. The resnltant 
simnlated TGFs have now nndergone the same processes as the real data. 

To generate a synthetic TGF incident on GBM with a photon flnence of n ph/cm^, 
the following steps are taken. First, a set of DRMs is created nsing the known sonrce 
location. Next, the appropriate annnlns from which to draw photons is selected. For each 
GBM detector, the nnmber of incident photons (nj) is given by the prodnct of n and the 
geometric area of the detector (Ag), the projection of the snrface area of the detector in 
the direction of the sonrce location. The desired nnmber of photons, nj, is randomly drawn 
from the sonrce distribntion and folded throngh the DRM for that detector. 

To fold a single photon with a DRM, the DRM must hrst be converted to a probability 
distribution by dividing by Ag. The appropriate channel column in this matrix is then be 
selected. The sum of this column is the probability that the photon deposits any energy 
in the detector. A random number is drawn and if it is less than the summed probability, 
a random channel is drawn from the column; otherwise, the photon is discarded. The 
background is simulated as a homogeneous Poisson process, with a rate set to the mean of 
the observed background rate distribution for the Nal and BGO detectors, which is 1.1 and 
1.6 kHz respectively. Each count is assigned a corresponding energy value which is randomly 
drawn from the observed background count spectrum distribution. The background is then 
combined with the source counts and passed through an instrumental hlter that includes 
the effects of both dead time and pulse pile-up and is based on Ghaplin et al. [2^ . 

As the off-axis distance is increased, the photon huence incident on GBM is normalised 
to match the observed counts. This is achieved by converting the observed count huence to 
photons/cm^ by dividing by the ehective area. For this step, we consider only the observed 
huence in the BGO detectors, as for these, the inhuence of dead time on the observed huence 
has been studied in detail (e.g. Tierney et al. |23], Briggs et al. [2H])- Nal detectors suher 
from a greater ehective dead time due to the higher proportion of overhow counts compared 
to the BGO detectors (overhow counts, which are above the maximum digitised energy, 
incur a 10/rs dead time). 

For each real TGF in our sample, a 15 km source model with a huence of 1 ph/cm^ was 
folded through the response of each BGO detector. The ehective area is given by the ratio 
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FIG. 1. Evolution of the duration (crosses) and hardness ratio (circles) of the RREA simulations 
at an altitude of 565 km, as a function of the TGE source to sub-satellite distance for a 15 km 
altitude wide beam source model. The duration is measured using the Tso, the time interval in 
which 50% of the flux occurs, starting and ending at 25% and 75% levels. The hardness ratio is 
given by the number of events with energy greater than 300 keV divided by the number of events 
with energy less than 300 keV. As the source off-axis distance is increased, a clear elongation in 
time and spectral softening is visible. This is due to the increased Compton scattering experienced 
by the photons as they propagate through a greater integrated density of atmosphere. 

of the observed counts from the DRM-folded model to the number of incident photons. The 
incident photon fluence of each TGF is then calculated by dividing the number of observed 
counts in the data by the effective area. Applying this method to the entire sample allows 
the approximate source fluence in ph/cm^ as a function of the TGF source to sub-satellite 
distance to be determined. This is shown in Figure]^ The fluence is binned into 50km bins 
to match the simulations and the average fluence in each bin found. In the absence of dead 
time, these values could be used as the incident photon fluence in the simulations. However, 
the dead time is non-negligible and a correction factor must be applied. 

The actual percentage loss depends on the incident flux, but losses of 60 % during the 
peak emission have been estimated for a very bright TGF |2H]. For simplicity, we assume an 
average loss of 40 %. The corrected fluence in each 50 km bin is used as the incident photon 
fluence in the corresponding simulation bin. 
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FIG. 2. Photon fluence in BGO detectors as a function of the source to satellite distance. These 
values have not been corrected for dead time losses. The horizontal grey line represents the cut 


that was applied to ensure robust statistics for the analysis (see § IV). 


The REAM code generates an initial source electron distribution that is created instan¬ 
taneously. In reality, the electron distribution will have a time dependence (e.g. Dwyer 
[29], Celestin et al. [30|). To add this feature to our simulations, we artificially smear the 
simulated photon arrival times at the spacecraft with a Gaussian distribution. As there 
are limited predictions in the literature for the timescale of the variation at the source, 
we consider 4 empirical smearing distributions, with standard deviations (a) of 25, 50, 75 
and 100 /rs respectively. These values were selected to be representative of the durations 
observed at spacecraft altitudes. 


IV. ANALYSIS 

In this work we analyse hundreds of real, and thousands of synthetic TGFs. It is desirable 
that we have an analysis method that is as objective and automatic as possible. The most 
problematic step is determining the point at which the TGF begins and ends. Once the 
interval is chosen, auxiliary information such as the hardness ratio can then be calculated 
over this interval. 
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FIG. 3. Results of the Bayesian Block Algorithm for TGF081123. The top panel shows the time 
prohle, with the integrated counts in 20 fis bins shown in grey. The solid line is the optimum 
representation of the data. The bottom panel shows the probability distribution for each change 
point. 


In previous analyses of GBM data, the Tgo measure has been used to define the time 
interval and duration of TGFs (e.g. Briggs et al. [28]). Adopted from the study of GRBs, 
this takes the duration of an event as the time interval in which the fluence rises from 
5% to 95% [31]. The time interval in which the counts are accumulated is dehned by the 
user, and thus requires human interaction. Due to the low number of counts detected, and 
general trend for the intensity of TGFs to fall with time, determining the times at which 
the TGF is indistinguishable from instrumental background can be difficult. Ref. [18], uses 
the related T 50 measure as it is less susceptible to uncertainties caused by low count rates 
and background. However, the T 50 is not necessarily representative of the TGF duration. 
For these reasons, we do not use the Tgo or T 50 measures in this work. 

To select the time intervals corresponding to the TGF for both actual observations and 
simulations, we employ the Bayesian Block Algorithm (BBA). This is a non-parametric 
algorithm that hnds the optimal segmentation of data [32] • It is designed to address the 
general problem of detecting and characterising variability in time series data and can be 
applied to both time tagged and binned data. The data are divided into discrete segments 
or blocks, in which the rate is modelled as a constant. In practice, this translates to binning 
the data into non-uniform bins of common rate. This algorithm is frequently used in X-ray 
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and gamma ray astrophysics (e.g. Buehler et al. [33]). A brief overview of the algorithm 
follows. 

The algorithm iterates over the data, adding in one data point with each iteration. As 
each data point is added, all possible segmentations of the data are tested. The segmentation 
which maximises the htness is chosen, the expression for which depends on the data type. 
For GBM TTE data, the relevant expression is given by Eq. 19 in Scargle et al. [32] . The 
nnmber of blocks is not explicitly set, bnt is inflnenced by dehning a prior distribntion 
for the nnmber of blocks. Ideally this prior shonld assign higher probability to a small 
nnmber of blocks. The geometric prior, given by P{N) = p'j^biocks^ jg ^ged. As A^biocks is 
not known in advance of the analysis, the contribntion of the prior to the htness is 

introdnced. Using simnlations of pnre-noise time series, is related to a false positive 

probability po- By adjusting po, the prior distribution is adjusted. A complete description 
of the algorithm can be found in Scargle et al. [32] . 

Using Po = 0.05, BBA was run for each TGF. The blocks corresponding to the TGF were 
selected by comparing them to the background rate. The time interval of these blocks was 
then taken as the duration of the TGF (Tbs). Figure shows an example of how the BBA 
technique is applied to a TGF. 

In comparison to Tgo, Tbb is more conservative while also being less subjective. To com¬ 
pare the results from the Bayesian Block analysis to the Tgo measure, we use the intersection 
of the sample used in this work and that used in Gonnaughton et al. [IH]. For 158 common 
TGFs, the Tgo values are plotted against the corresponding Tbb in Figure (although only 
T 50 values were published in Gonnaughton et al. ra. Tgo values were also produced by the 
same analysis). To quantify the degree of correlation between the two measures, we use 
the Pearson product-moment correlation test, obtaining a coefficient of 0.47 which implies a 
moderate degree of linear correlation and indicates that the two measurements are broadly 
consistent. 

For each TGF, the duration was taken as the time interval dehned by Tbb- The hardness 
ratio (HR) of an event is dehned as the ratio of counts above an energy threshold to those 
below it. This measure is useful as it can be used as an identiher for spectral evolution. 
However, the HR is heavily instrument dependent as it is based on count data, without 
any deconvolution to incident photons. Following previous studies (e.g. Grefenstette et al. 
0), we have adopted the convention of dehning soft counts as those with measured energies 
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FIG. 4. Comparison between duration calculated using Bayesian Blocks (Tbb) and Tgo values 
for the intersection of TGFs used in this work and Connaughton et al. |18j . The shaded region 
indicates the 68% containment region. 



<300keV. For each TGF, the HR (Fl> 3 oo/.F< 3 oo) was calculated over the time range Tbb- 

Following Grefenstette et al. [8], we characterise the spectral evolution of a TGF by the 
delay between the counts above and below 300 keV. This is calculated over the interval 
Tbb by first finding the average arrival time of each component and then calculating the 
difference. The error on each component was taken as the standard error on the mean. The 
error on the delay is the quadrature sum of each. 

As the source distance is increased there is a decrease in the observed fluence. To mitigate 
the effect of this on the analysis and to ensure robust statistics, TGFs with an observed 
fluence less than 0.08 ph cm“^ and source distances greater than 500 km were discarded. 
The fluence cut is shown as a horizontal grey line in Figure 

The simulated TGFs were analysed in an identical fashion to the data. To facilitate the 
comparison of this analysis to the data, the distributions of derived parameters (duration, 
delay, HR) were analysed and fit with a Gaussian for each set of 1000 simulations. The 
mean and standard deviation of each fit could then be compared to the data directly. 
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V. RESULTS 


The distribution of observed delays in TGFs is shown in Figure The observed values 
range from -20 to 80 fis with a mean of ~ 24 /rs. It is tempting to compare this to the value of 
28 fis obtained in the RHESSI analysis of 431 stacked TGFs [H], however care must be taken 
as as the delays are calculated from counts detected in various energy bands, measurements 
that are detector dependent. The RHESSI data were also corrected for dead time losses. 

The delay as a function of the TGF source to sub-satellite location distance is shown in 
Figure for the data and simulations. The data exhibit no significant variation with offset. 
The simulations which have undergone no smearing show no significant variation and are 
consistently longer than the mean of the data. The four electron timescales considered for 
the smeared simulations are all broadly consistent, and up to 200 km, are consistent with 
the mean of the data. After this, similar to the zero smearing cases, the simulated delays 
are consistently longer than the data. 

The hardness ratio as a function of the TGF source to sub-satellite location distance is 
shown in Figure for the data and simulations. A clear trend of increased softening with 
offset is visible in the data. The simulations which underwent no smearing are completely 
inconsistent with the data. Simulations with a 25 /rs smearing factor are inconsistent with 
the data up to 150 km. After this point they follow the trend but are significantly softer than 
the data mean. The simulations with longer smearings (50, 75 and 100 /is) are consistent 
with each other, and agree best with the data, but are also systematically softer than the 
mean of the data. Possible explanations for this discrepancy are discussed in the following 
section. 

The duration as a function of the TGF source to sub-satellite location distance is shown 
in Figure for the data and simulations. The data exhibits a decreasing duration with 
increasing source distance, the opposite of the predicted relation. This is likely a consequence 
of the decrease in fluence with increasing source distance (Figure]^, as it will be more difficult 
to distinguish the TGF from the background, resulting in a shorter observed duration. The 
durations from the simulations exhibit little variation with increased offset, and encompass 
the range of observed durations within their spread. The results from the 50 /rs smearing 
agree with the mean of the data up to 300 km, after which on they are consistently longer 
than the data. 
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FIG. 5. Distribution of the observed delay between soft and hard counts. 

VI. DISCUSSION 

Assuming a wide beam model and a single altitude source, an increase in the off-axis 
distance is expected to cause a corresponding increase in duration and decrease in spectral 
hardness (see Figure [^. This temporal elongation is not observed in either the data or the 
simulations, which is likely to be due to the decrease in the incident fluence as the off-axis 
distance is increased. In contrast, the predicted softening with increased off-axis distance is 
evident in the data. The simulations with smearings less than 50 /rs are generally inconsistent 
with the mean of the data, implying that the source electron distribution for the majority 
of TGFs observed by GBM is not necessarily created instantaneously, and likely varies on 
timescales > 50 /rs. 

It has been hypothesised that short TGFs observed by GBM (T 50 < 50 /is) are the result 
of instantaneous sources, and those with longer durations are the result of a superposition of 
these pulses (TB] . This analysis would seem to challenge the former, as the spectrum obtained 
from a short electron pulse (< 25/is) is inconsistent with the data. This implies that the 
temporal elongation due to Gompton scattering is not sufficient, and that an intrinsic time 
variation of at least tens of /rs is required at the source. 

Even for the largest smearing considered (100/is), the data exhibit harder spectra than 
can be explained by the simulations. We posit three possible explanations for this discrep¬ 
ancy in order of likelihood. The assumption of a single source altitude of 15 km may be 
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□ Data Mean ^ ct=0 B ct= 100 ^xs 



FIG. 6. Delay between soft and hard counts as a function of the distance between the source 
and the Fermi sub-satellite position. The individual TGFs are plotted in grey, and the average 
in 50 km bins is plotted as cyan squares with the standard deviation plotted as an error bar. The 
values obtained in the simulations are plotted in orange. For the sake of clarity, the simulated 
values are offset from the data in 5 km steps. To indicate the spread of the simulated values, the 
standard deviation of the fit is plotted as an error bar. The different smearings for the simulations 
are indicated by the marker, circles for no smearing and squares for 100 /xs. 

incorrect. If instead there are a range of source altitudes, then the evolution of HR with off- 
axis distance would be considerably broadened. TGFs with harder spectra may be simply 
explained by a higher source altitude, as the photons will undergo less Compton scattering 
due to the reduced integrated atmosphere traversed. Alternatively, the RREA model used 
in this work may be too basic, and more complex models which include an evolving elec¬ 
tron source distribution (e.g. Dwyer [21]) may be required to explain the range of observed 
values. Finally, ‘tilted’ beams may be invoked. In this scenario, the alignment of the E 
held at the source is not directly parallel to the vertical, with the result that the beam of 
gamma rays is correspondingly tilted. In such a case, the off-axis HR would be harder than 
expected due to the greater proportion of high energy photons. However, the same general 
trend of decreasing hardness with increased off-axis distance would remain, as the tilted 
beams would traverse a greater depth of atmosphere than the vertical equivalent. 

It also possible that the explanation is a combination of the hrst two points above. TGFs 
close to the satellite nadir could have a range of source altitudes with primarily vertical 


15 


























FIG. 7. Hardness Ratio as a function of the distance between the source and the Fermi sub¬ 
satellite position. The data are plotted in the same fashion as Figure The expected softening 
with increased offset is evident. The values derived from the simulations with no smearing are 
significantly softer than the data. Simulations with longer smearing timescales (>50//s) are more 
representative of the data, but still tend to underestimate the hardness. In the interests of clarity, 
only the results from the 100 //s simulations are overplotted. 


□ Data Mean ^ ct=Q ^is ^ ct= 50 ^iS B cT=100^iS 



FIG. 8. Calculated duration as a function of the distance between the source and the Fermi 
sub-satellite position. The data are plotted in the same fashion as Figure The durations from 
the simulations exhibit little variation with increased offset, and encompass the range of observed 
durations within their spread. The results from the 50 ^s smearing agree with the mean of the 
data up to 300 km. 
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beams. As the source distance is increased, TGFs with higher altitude sources and/or 
tilted beams would be preferentially detected. Accounting for this possible selection effect 
would introduce degeneracies that cannot be mitigated with our sample size. They will be 
investigated in future analyses using simulations and a larger sample size. 

Our results imply that the majority of TGFs observed by GBM must have source electron 
distributions that vary on timescales of at least tens of /rs. Based on Figure we suggest 
that 100 /is is close to the upper limit of the source electron distribution variation timescale, 
with a value of 50 /is being a likely mean. Of the hve empirical smearing distributions 
considered, only the 50 /is Gaussian is consistent with the observed temporal and spectral 
properties of the data (i.e. delay, duration and hardness). This consistency extends up 
to 200 - 250 km, at which point the simulated times become longer and the spectra become 
softer. This could be attributed to insufficient statistics beyond this source-satellite distance 
(see Figure]^. 

In this work, we have performed a comprehensive study of RREA simulations in compar¬ 
ison to GBM observations. The observations exhibit a clear softening with increased source 
distance, in qualitative agreement with theoretical predictions. However, a quantitative 
analysis shows that the observed spectra can be harder than that predicted by the basic sin¬ 
gle source-altitude model. Simulations performed with an electron source timescale of 50 /rs 
are most consistent with the temporal and spectral properties observed in the data. We pro¬ 
pose that the source electron distributions of TGFs observed by GBM vary on timescales 
of at least tens of microseconds, with an upper limit of ?«100 /rs. Drawing more concrete 
conclusions is limited by the low number of TGFs at larger off-axis source distances and 
the simplifying assumptions of the simulations. These assumptions, which are common in 
the literature, include a single hxed source altitude of 15 km and a vertically oriented E 
held. The effect of varying the source altitude and beam orientation will be investigated in 
a future work. 
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